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Abstract. We report on the next-to-leading order QCD calculation for e + e — — > 
4 jets. We explain some modern techniques which have been used to calculate 
the one-loop amplitudes efficiently. We further report on the general purpose 
numerical program "Mercutio", which can be used to calculate any infrared safe 
four-jet quantity in electron-positron annihilation at next-to-leading order. 



1. Motivation: LEP physics 

QCD four-jet production in e + e~ annihilation can be measured at LEP and can be 
studied in its own right. First of all, e + e~ — > 4 jets is the lowest order process 
which contains the non-abelian three-gluon-vertex at tree level and thus allows for an 
measurement of the colour factors Cf, Ca and Tr of QCD. This in turn may be used 
to put exclusion limits on light gluinos. Furthermore the QCD process is a background 
to W-pair production, when both W's decay hadronically and to certain channels of 
the search of the Higgs boson like e + e~ — > Z* — > ZH — > 4 jets. The one-loop matrix 
elements required for an NLO study of four-jet production are also an essential input 
for an NNLO calculation of three-jet production. The latter one would be needed to 
reduce theoretical uncertainties in the extraction of the strong coupling at the Z-pole. 

In general leading-order calculations in QCD give a rough description of the 
process under consideration, but they suffer from large uncertainties. The arbitrary 
choice of the renormalization scale gives rise to an ambiguity, which is reduced only 
in an next-to-leading order calculation. Furthermore the internal structure of a jet 
and the sensitivity to the merging procedure of the jet algorithm are modelled only 
in an NLO analysis. Both uncertainties are related to the appearance of logarithms, 
ultraviolet in nature in the first case, infrared in the latter, which are calculated 
explicitly only in an NLO calculation. 

A NLO calculation proceeds in two steps: First, one needs the relevant 
amplitudes, in our case e + e~ — > 4 partons and e + e~ — > 5 partons at tree level and 
e + e~ — > 4 partons at one loop. Among these, the one-loop amplitudes are the most 
complicated ones and we will comment on their calculation in the next section. The 
second step requires setting up a numerical Monte Carlo program which has to deal 
with infrared divergences. We will focus on this point in the third section. In the last 
section we will give some numerical results. 
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2. One-loop amplitudes 



We used a variety of modern techniques in order to calculate the one-loop amplitudes 
efficiently. These include colour decomposition, where amplitudes are decomposed 
into simpler gauge-invariant partial amplitudes with definite colour structure, and 
the spinor helicity method, which consists in expressing all Lorentz four-vectors and 
Dirac spinors in terms of massless two-component Weyl-spinors. Their use divides 
the task into smaller, more manageable pieces. Also a decomposition inspired by 
supersymmetry proved to be useful, where the particles running around the loop 
are reexpressed in terms of supermultiplets. In a second step the cut technique and 
factorization in collinear limits are used to constrain the analytic form of the partial 
amplitudes. 

As an example we explain in more detail the cut technique which is based 
on unitarity. To obtain the coefficients of the basic box, triangle or bubble integrals 
one considers the cuts in all possible channels. Each phase-space integral is rewritten 
with the help of the Cutkosky rules as the imaginary part of a loop amplitude. The 
power of this method lies within the fact, that on each side of the cut one has a 
full tree amplitude and not just a single Feynman diagram. This method allows one 
to reconstruct the one-loop amplitude up to terms without an imaginary part. The 
remaining terms were obtained by examining the collinear limits. 

For the reduction of tensor pentagon integrals we used a new reduction algorithm 
, based on the Schouten identity and Weyl spinors, which does not introduce artifical 
Gram determinants in the denominator. 

The one-loop amplitudes for the first subprocess e + e~ — > qqQQ were calculated 
in refs. J3], ^| and the amplitudes for the second subprocess e + e~ — > qqgg in refs. 
P, a. The calculations of the two groups agree with each other. 



3. Numerical implementation: Mercutio 

The second major part of a general purpose NLO program for four jets is coding the 
one-loop amplitudes and the five parton tree-level amplitudes in a numerical Monte 
Carlo program. At leading order the task is relatively simple: One parton corresponds 
to one jet. At NLO however, a jet can be modeled by two partons. At NLO the cross 
section receives contributions from the virtual corrections and the real emission part. 
Only the sum of them is infrared finite, whereas when taken separately, each part gives 
a divergent contribution. Several methods to handle this problem exist, such as the 
phase-space slicing method 0, the subtraction method || and the dipole formalism 
||. We have chosen the dipole formalism. Within the dipole formalism one subtracts 
and adds again a suitably chosen term: 




The approximation term do~ A has to fullfill the following two requirements: First, da A 
must be a proper approximation to do~ R , with the same pointlike singular behaviour in 
D dimensions as da R . Secondly, da A must be analytically integrable in D dimensions 
over the one-parton subspace leading to the soft and collinear divergences. 

Let me now turn to the details of the Monte Carlo integration. The heart of 
any Monte Carlo integration is the random number generator. Among other things, 
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it should have a long period and should not introduce artifical correlations. As the 
default random number generator we use 

s 4 = (si-24 + Si_ 55 ) mod 2 32 . (2) 

It was proposed by Mitchell and Moore and has a period of 2^(2 55 — 1), where < 
/ < 32. Massless fourmomenta are generated with the help of the RAMBO-algorithm 
[jTo|| . This algorithm generates events with a uniform weight. Adaptive importance 
sampling is implemented using the VEGAS-algorithm @. A naive implementation 
of the dipole formalism will give large statistical errors when performing a Monte 
Carlo integration over the real corrections with dipole factors subtracted. In order to 
improve the efficiency of the Monte Carlo integration we remap the phase space to 
make the integrand more flat. A simplified model for the term do~ R — do~ A would be 

F = hjM-m ( 3) 





where /(0) = g(0) is assumed. f(x)/x corresponds to the original real emission part 
with a soft or collinear singularity at x = 0, g(x)/x corresponds to the subtraction 
term of the dipole formalism. Eq. (J3j) can be rewritten as 

9m,„ 

F= J dx f{x) ~ 9(x) + J dy(f( e y)-g(ey)), (4) 

In J/min 

where y m i n is an artificial parameter separating a numerically dangerous region from 
a stable region. Using the Taylor expansion for f(x) — g(x), one sees that the first 
term gives a contribution of order 0(y m i n ). In the second term the l/x behaviour has 
been absorbed into the integral measure by a change of variables y — hix, and the 
integrand tends to be more flat. It should be noted that there is no approximation 
involved. 



4. Numerical results 

Numerical programs for e + e~ — > 4 jets have been provided by four groups: MENLO 
PARC U), DEBRECEN |J), EERAD2 @ and MERCUTIO |§J|. Various cross- 
checks have been performed among these programs and they agree within statistical 
errors. Here we report on the numerical program "Mercutio" , which was written in 
C++. The four-jet fraction is defined as 

Ri = ?±±L. (5) 

O~tot 

The values obtained for the four-jet fraction for the DURHAM algorithm with 
Vcut — 0.01 for various energies are given in table ^. The decrease with energy is 
mainly due to the running of the strong coupling. 

With the numerical program for e + e _ —>■ 4 jets one may also study the internal 
structure of three-jets events. One example is the jet broadening variable defined as 

Bjet = ^5?fm (6) 

J jets ' 
J a 

Here is the momentum of particle a transverse to the jet axis of jet J, and the sum 
over a extends over all particles in the jet J. The jet broadening variable is calculated 
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Vq 2 


R\° 


ftNLO 


135 GeV 
161 GeV 
172 GeV 
183 GeV 
189 GeV 


(2.98 ±0.01) • 10"2 
(2.65 ±0.01) • 10- 2 
(2.53 ±0.01) • 10- 2 
(2.48 ±0.01) • 10- 2 
(2.44 ±0.01) • 10- 2 
(2.42 ±0.01) • 10" 2 


(4.72 ±0.01) • 10" 2 
(4.12 ±0.01) ■ 10- 2 
(3.89 ±0.01) ■ 10- 2 
(3.81 ±0.01) ■ 10- 2 
(3.73 ±0.01) ■ 10- 2 
(3.69 ±0.01) • 10" 2 



Table 1. The four-jet fraction at LO and NLO for the DURHAM algorithm 
with y cu t = 0.01 and various energies. 




Figure 1. The Bj et distribution at NLO (diamonds) and LO (crosses). 



for three-jet events denned by the DURHAM algorithm and y cu t =0.1. This choice 
is motivated by a recent analysis of the Aleph collaboration [|l6| . Figure [l] shows the 
distribution of the jet broadening variable. 
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